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We study magnetic fluctuations in a system of interacting spins on a lattice at high temperatures 
and in the presence of a spatially varying magnetic field. Starting from a microscopic Hamiltonian 
we derive effective equations of motion for the spins and solve these equations self-consistently. We 
find that the spin fluctuations can be described by an effective diffusion equation with a diffusion 
coefficient which strongly depends on the ratio of the magnetic field gradient to the strength of 
spin-spin interactions. We also extend our studies to account for external noise and find that the 
relaxation times and the diffusion coefficient are mutually dependent. 

PACS numbers: 



I. INTRODUCTION 

Recent advances in magnetic imaging techniques, as 
well as the development of novel types of electronic de- 
vices that utilize electronic spin (rather than charge) as 
an information carrier, have renewed interest in under- 
standing mechanisms of spin noise and spin relaxation. 
While conventional experimental methods, such as nu- 
clear or electron spin resonance and related techniques^, 
probe the temporal evolution of spin correlations, they 
typically do not provide much information on spatial cor- 
relations between neighboring spins. On the contrary, 
the new approaches to spin resonance, such as magnetic 
resonance force microscopy (MRFM), combine capabili- 
ties of the usual magnetic resonance techniques with the 
sensitivity of atomic force microscopy. That is, one can 
now observe not only the time (frequency) dependence 
of spin correlations, but also their spatial dispersion with 
an atomic-scale resolution. Hence, there is a clear need 
to develop theoretical tools for the description of such 
correlations in systems of interest, that is, in systems of 
interacting spins. 

The spatial correlations in interacting spin systems are 
believed to be controlled by the so-called flip-flop pro- 
cesses. That is, two neighboring interacting spins can 
exchange magnetization, i.e., the values of their spin com- 
ponents can change by ±1/2, so that the total spin of 
the pair is conserved. Such exchange gives rise to the 
diffusion of spin magnetization, provided the dynamics 
of the flip-flops is Poissonian^. Typical calculations of 
the effective diffusion constant utilize the method of mo- 
ments, where the line- width is approximated by a gaus- 
sian or lorentzian shaped. Such approximations are not 
very well controlled. More recently several types of clus- 
ter/cummulant expansions have been proposed in con- 



nection with the problem of decoherence of localized 
electronic spins caused by the fluctuations of nuclear 
spins^— . In that problem though, the decoherence of 
electronic spins occurs on a timescale small compared to 
the typical nuclear timescale, which justifies the use of 
cluster expansions in the description of fluctuations in 
the nuclear subsystem. 

In this paper we study correlations between spatially 
separated spins in the opposite, long time regime. Such 
a regime is specifically relevant to the MRFM technique, 
which utilizes (micro) mechanical cantilevers with ferro- 
magnetic tips to probe magnetic fluctuations in the un- 
derlying samples. We propose an approach based on the 
Markov approximation, similar to the frequently used 
Bloch-Rcdficld approximation^ in the theory of open 
quantum systems. That is, we consider all possible pairs 
of interacting spins, while other spins ^ are 
treated as an environment, providing finite line- width 
for the flip-flop transitions through fluctuating magnetic 
fields (see Fig. [1] for a cartoon visualisation of these ap- 
proximations) . A self-consistency is then established be- 
tween the flip-flop rates and the line-width so that our 
approach can be viewed as a sort of dynamical mean field 
approximation. We argue that our method is well justi- 
fied, in particular, in the presence of an external strongly 
non-uniform magnetic field, which introduces separation 
between the timescales of the flip-flop rates and the corre- 
lation time for the fluctuations of the effective magnetic 
fields. Note that such non-uniform magnetic fields are 
intrinsic to the MRFM setups, where field gradients are 
used to address specific spins located within the so-called 
resonance layer. 

Our paper is organized as follows. In Section|TI]we de- 
scribe a general formalism that can be utilized to study 
spin-spin correlations for a broad class of spin Hamilto- 



2 



Hamiltonian 



(b) 




FIG. 1: (a) shows a collection of spin-half particles on a rigid 
lattice in a nonuniform external magnetic field. The quantity 
of interest in this work is the rate at which spin flip-flops 
occur. Our model is displayed pictorially in (b), where the 
neighbouring sites of i and j are replaced by a fluctuating 
bath. The flip-flop rate is then calculated such that it is 
consistent with the fluctuations of the bath. 



nians, e.g. Eq. ([T]). We derive effective equations of 
motion for the magnetization, e.g. Eq. (|12[) . which has 
the form of a stochastic master equation. In doing so wc 
use methods developed in connection with studies of dif- 
fusion in classical lattice gas models^ as well as in the 
theory of open quantum systems^. The equation of mo- 
tion is supplemented by a self-consistency equation, Eq. 
(|T3"|) , which relates the rates in the master equation to the 
correlation function evaluated from the master equation 
in terms of the rates. In Section Hill we look specifically 
at the Heisenbcrg model on a cubic lattice in the pres- 
ence of a spatially non-uniform external magnetic field. 
We find that the flip-flop rates are strongly suppressed by 
the field gradient in the limit when the field gradient sig- 
nificantly exceeds the spin-spin interaction constant. In 
Section ITVl we study the influence of spin-relaxation pro- 
cesses on spin flip-flops and derive the effective master 
equation for the magnetization in the presence of exter- 
nal noise sources acting on the spins. Our main result 
of that section is that, while the field gradient suppresses 
the flip-flops, the noise may actually enhance these rates; 
see Eq. ([53)) and corresponding discussion. Finally, in 
Section [V] we discuss the validity of our approximations 
and summarize the results. 



II. MODEL AND GENERAL SOLUTION 

We consider a system of spin-half particles on a lattice, 
interacting with each other according to the following 
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(1) 

where a k = (crjj: ± icr\)/2, k = (k x ,k y ,k z ), and <r£ are 
Pauli matrices, a = x,y 1 z. The index i in the first sum 
runs over all lattice sites, while the notation in the 
second sum indicates the summation over all pairs of lat- 
tice sites. The external magnetic field Bi is assumed to 
be non-uniform in space. The spin-spin interaction is 
isotropic when jjL = J^.. The equation of motion for cr| 
is 

id t at = [af,H] = 4 £ 4(<7+<r fc " - (2) 

k^i 

In Eq. @ and in the following we set H = 1. Next 
we consider the equation of motion for <jfa~^. After a 
straightforward calculation we obtain 



id t {aial)=2 5B^a+al + Ji{at-al) (3) 

" Jik v 



n^{i,k} 

where SB^f is the difference between effective magnetic 
fields at sites k and i, 
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This difference consists of a constant part; 

5Bki = Bk — Bi, 



(4) 



(5) 



and a part which fluctuates (due to spin flips at nearby 
lattice sites); 



sB R k r(t)= J2 [ J L< 



Til Z 

'J ni n 



(G) 



The larger the number of individual spins contributing 
to 8Bfa Ct , the more rapidly fluctuating this quantity be- 
comes. Hence, for systems with sufficiently long-range 
interactions or high dimensionality 8B^f ct fluctuates very 
rapidly. 

From Eq. © we see that the expectation value of 
<7^~<7^T contains a prefactor 

A^MHe 2 ^'^'), (7) 

related to the Larmor precession of spins around the ef- 
fective magnetic field at sites i and fc. The fluctuating 
component of the effective-magnetic-field [see Eq. ©] 
causes the precession frequencies at each site to vary. 
Moreover, if the effective magnetic fields at sites i and 
k are large, and the number of spins contributing to the 
fluctuating component of the field [see Eq. ^] is much 
greater than one, then from Eq. [or more specifically, 
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the prefactor shown in Eq. (J7J], we would expect the Lar- 
mor precession frequency of cr^cr^ to be very fast (com- 
pared to the dynamics of the individual <rf operators) 
and fluctuate rapidly. Following this logic, we see that 
the summation of terms Pn a k ano - a t G n m Eq. © is 
essentially a summation over a rapidly fluctuating ob- 
ject, and will statistically self-average to zero (provided 
a sufficiently large number of spins contribute to <5Bj^ uct ). 

A similar approximation is very common in the theory 
of open quantum systems, where it is known as the sec- 
ular or Bloch-Redfield approximation^.. As in the case of 
open quantum systems it relies on the assumption that 
the off-diagonal elements of a system's density matrix p 
are small cither due to large splittings between the adja- 
cent energy levels or due to rapid fluctuations from the 
heat bath. In the present case the fields SB^ ct (t), play 
the role of the heat bath operators and must treated self- 
consistently, to which we now focus our attention. 

By integrating Eq. ([3]) (with the summation on the 
right-hand-side neglected) we obtain 



o+o-fc (t) ~ -Ui \ dsA ik (t,s)[af(s) - <rf( S )] (8) 
Jo 

+A ik (t, 0)cj fc , 

where the last term is due to the initial condition of the 
operator ci k = <jfa^{t = 0). In the high temperature 
limit the system is disordered and therefore it is natural 
to assume that the expectation value of crj cfZ is random, 
with ((o+o£» = and ({vf a* a$ a„)) = (1/4)5^8^, 
provided i ^ k and i' ^ k' . Here the double bracket 
stands for averaging over the ensemble of density ma- 
trices of the system as well as over a particular real- 
ization of the density matrix (set by a particular choice 



much faster than the evolution of [of (s) — <7|(s)], we can 
employ the Markov approximation, and set s = t which 
removes the latter term from the integral in Eq. © to 
give 



-iJi 



ik , ds(A ik (t,s))[af(t)-am] (H) 
o 



+A lk (t,0)c ik . 



Now that we have a formal solution for <jt a k W it is pru- 
dent to substitute the expression back into the sum in 
Eq. © which was originally ignored in deriving Eq. (jlip . 
In doing so we wish to find an inequality which quantita- 
tively ensures the summation term is small compared to 
all other terms in Eq. ([3]). The details of this calculation 
are straightforward (sec Scction|V]for further discussion) 
and one finds Ji k <C Ti k (where is the rate at which 
flip-flops occur and is calculated below) is a sufficient 
condition to ensure the summation in Eq. (|3]) remains 
small. 

We now substitute Eq. (fTT|) into Eq. @, to give 
dt(am) =E!> [<°i(*)> - KW)] +&(*)• (12) 

The averages in Eq. ([12]) are taken with respect to a par- 
ticular realization of the systems density matrix, but not 
over the ensemble of the density matrices. The coeffi- 
cient, Tj-fc, represents a rate at which spin flip-flops occur 
between sites j and k (these can only occur when sites j 
and k have opposite spin). The expression for this rate 
is given by 



of the initial condition), i.e., u k 

(Wt°k)) = (Tr(a+a,Tp)) P , etc. 

We wish to substitute Eq. ((8]) into Eq. ([2]) to obtain a 
closed form equation for o\(t). This can be significantly 
simplified if we replace the rapidly fluctuating quantity, 
Ai k (i, s), in the integrand in Eq. © by its average value. 
This approximation is in a perfect agreement with our 
assumption regarding the separation between time scales 
for the dynamics of the local fluctuating magnetic field 
at site i, and components of the individual spin at site i. 
We make the assumption that, by virtue of the central 
limit theorem, the random variable SBfj? is Gaussian; 

(A ik (t,s)) = e HBi-B k )(t-s) e -2^ K^in-^dnd^ 

= e 2i(B z -B k )(t-s) e -Af^- sl K ik (ji)(\t-8\-n)dn 

(9) 



i particular choice r°° 
= Tr(afa7p) and = 4 ( J jfc) / CX P 

J — oo 



2iSBj k s 



\s\ 

4/ K kj (ji){\a\-fi)dn 
lo 



ds, 



(13) 



where 



Kikin - r 2 ) 



(10) 



where we have used the quickly-decaying property of 
(Ai k (t, s)) to extend the upper and lower limits of the 
integral to ±oo. The final term in Eq. (fl~2]) represents 
the uncertainty with respect to the choice of the initial 
conditions of the system, and is given by 

&(f) =4z£ Jf k [c jk A kj (0, t) - c kj A jk (0, t)] . (14) 

Averaging over £i(t) corresponds to averaging over an 
ensemble of different density matrices (each density ma- 
trix being distinguished by a unique initial condition). 
Noting that (A ifc (0 > t)A w (0,f)) = (A w (f ,t)), and since 
Ai k {t' , t) is a rapidly fluctuating function of t — t', we can 
make the approximation; 



is the autocorrelation function of the fluctuating compo- 
nent of the magnetic field gradient between sites i and k. 
Moreover, since Eq. © (as a function of \t — s\) decays 



= 26 (t ~ ~ r 3k + s jk r ^ 



(15) 
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Together, Eqs. (|T2|) and (fT5|) obviously describe Pois- 
sonian dynamics of a coupled two-state system. Indeed, 
we could have obtained the same result if we had pos- 
tulated that the dynamics of a given spin (say, at site 
i) is controlled by its flipping rates — Ylk^ ikCri (l ~ <Tfc ) 
and I]fcf fei crfe(l - <Ji), where f ik = T ik + rj ik , with T ik 
and r\i k being the constant and fluctuating parts of the 
rate respectively. In this case £, = Y^ki^ik — Vki)-> c -f- 
Eq. (|14[) . Note that one can derive Eq. ([13]) for the rates 
Fjfe within a straightforward perturbative calculation, as 
shown in Appendix [X] There, we calculate the probabil- 
ity of a flip-flop for a pair of spins in the presence of an 
external fluctuating field (along the z-direction) . In the 
current section, we have simply assumed that this fluctu- 
ating external field has been created by the neighbouring 
spins coupled to this pair (see Appendix [A] for details). 

Equations ([12]) and (|15[) constitute a closed system of 
equations, which allows one to evaluate the correlation 
functions (af(t)a k (t')). For an arbitrary choice of spin- 
spin interaction constants jf k and and external fields 
Bi, the rates Ti k in Eqs. (TT2]) and (TT5) . though formally 
unkown, are expressed in terms of these correlation func- 
tions (o-?(t)a%(l/)) [see Eqs. (13), (ID), and ©]. By 
evaluating these correlation functions in terms of To-, 
one obtains a closed set of equations which one must 
solve self consistently for r^. This provides a way of 
solving for both the rates, Ti k and the correlation func- 
tions, (erf (t)crf(t')) for an arbitrary choice of interaction 

constants; jJL, and external fields; Bi. 

In the next section we will evaluate the and 
(erf (i)crf(f/)) for a simple choice of coupling constants 
given by the three dimensional, cubic, Heisenberg model 
with nearest-neighbor interactions. 

Before proceeding to this task we note that in the limit 
of large field gradient \B{ — B k \ ^> jf k , the integrand of 
Eq. ([T3]) rapidly oscillates and therefore the value of the 
integral decreases with the growth of \Bi — B k \. In the 
limit of vanishing rate (erf {t)a k (t')} ~ 6i k , we find 



K-ik{t - f) ~ K tk 



E (4 



Evaluating then, the Gaussian integral in Eq. (flU) we 
obtain 



ik 



cxp 



(16) 



Thus we predict the rate at which flip-flops occur, and 
therefore the rate at which spin diffusion occurs, is very 
small for \B t - B k \ > j} k . 

III. EXAMPLE: HEISENBERG MODEL 

We now consider a particular example; the Heisenberg 
model on a cubic lattice with an external spatially vary- 
ing magnetic field. The Hamiltonian of the system can 



be cast in the form 

H = J2 BfcK + J E < <+e„ • ( 17 ) 



where i = (i x ,i y ,i z ), r 4 = z^ax + iyay -\-i z az (a being the 
lattice spacing) , v = 1 , . . . , 6 enumerates the unit vectors 
which point to the nearest neighbors: e\(2) = ±x, e 3 (4) = 
±y and e 5 ( 6 ) = ±z, and finally a = x,y, z. We also 
assume that the external field varies linearly in space, 
B(r) = bov g where g is a unit vector which points in the 
direction of variation. The Hamiltonian (|17[) obviously 
belongs to the class of Hamiltonians defined in Eq. (1). 

The equation of motion for erf. is given by Eq. (TT2") . 
which, for the Hamiltonian in Eq. (T7) reads 

9 t «) = E F e, [«+eJ - (<>] + CrM (18) 



and the noise £ ri (t) is correlated according to Eq. () L5[) , 
which becomes 

(0) = 2 S(t-t')J2r ev (5 rirj -s Ti+BvTi ). (19) 



Eqs. (15) and (TTD) can be readily diagonalizcd by a 
Fourier transform method. Writing 



a 2 (k,w)e lr - k+ ^*, (20) 



2ir J_ w/a (2tt) 3 



where the k-intcgral is taken over the first Brillouin zone, 
(a cube with an edge 2tt/q), we obtain from Eq. (fT5|) that 



<|a*(k,c)| 2 ) = 



(|f(k,w)| s 



^ 2 + {E,r e „ [l-cos(ae„k)]} 



f, (21) 



with i = x,y,z and £(k, cj) being the Fourier transform 
of £ r ,(t), defined similarly to Eq. (3D). From Eq. (ID) 



(ll(k-w) 



2 r e „ [1 — cos (a e„k)] 



(22) 



and taking the inverse Fourier transform of Eq. ([22) . we 
obtain 

K(t)crf i ,(0))=e-*^ r -x 

i n .(2r ei t)in,(2r ea t)i n ,(2re 5 t) ) 

(23) 

where /n(,z) is the modified Besscl function of complex 
argument^ and n x = \i x — i' x \, etc. At sufficiently large 
distances (and times) Eq. (|23p describes (anisotropic) 
diffusion with diffusion constants D vv ~ r e „et 2 . 

The rates T ev arc yet to be determined. They can 
be found from Eq. (TT3") . Note that while for arbitrary 
direction of the field gradient g the rates T Gv , T e , differ 
from each other, they are equal (T eu = T) for g = g = 
(l/v / 3)(x + y + z), i.e., when the field gradient points 
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FIG. 2: Numerical solution of the integral Eq. (|24|1 . showing 
the rate Y as a function of magnetic field gradient bo- 



along the main diagonal of the cube formed by the unit 
vectors x, y and z. In this case Eq. (fT5j) reduces to 



T = 4J 



where K(fi) is the correlation function of Eq. (JTUJ), which 
is now independent of the indices i and k, due to our 
convenient choice of magnetic field gradient direction g, 
which makes the diffusion process isotropic. K([i) can be 
easily expressed in terms of (of. (m) "^ (0)): 

k{h)=\j 2 (eEk>)<,(°)> 



,(0)) 



(25) 



where we have chosen to calculate K between sites r 2 = 
(0, 0, 0) and r 3 = e x (and then relied on the isoptropy of 
all directions in the lattice). Using Eq. (|2"3"|) we obtain 
K(p) = ±J 2 e~ 6I >/(2i», where 

f(x) =5tf(x) + Wh{x)ll{x) + Ul{x)I 2 (x) 

- Ul{x)h{x) - 8lf(x) - I2h{x)h{x)h{x) 
-Il{x)h{x). (26) 

Substituting this new found expression for K(fi) into Eq. 
(|2~4"|) we obtain an integral equation for T. One can solve 
this integral equation numerically to find T / J as a func- 
tion of bo/ J (see Appendix [B]), the results are shown in 
Fig. [5] For bo/ J 3> 1 the value of T is consistent with 
Eq. (|T^|) . which for the present case reduces to 

r ~4J^ e - h o/( 6a/2 ). (27) 

The analytic solution is also shown in Fig [2] for compar- 
ison. We find the analytic and numerical solutions are 
equal beyond bo > 10J. 



IV. INFLUENCE OF RELAXATION 
PROCESSES 

In this section we consider the influence of external 
noise on the spin-spin correlation function. We consider 
a model described by the Hamiltonian 



H = H 



(28) 



where H is given by Eq. ([TJ and r)f(t) is a fluctuat- 
ing magnetic field. The index i runs over lattice sites, 
and a = x,y, z. In reality such a field may arise due to 
phonons (for instance in semiconductors) or conduction 
electrons (for instance in metals). We will assume that 

(V? (.tfVj (*')) = <W 8ij A(f — t'), where A(t) is some even 
function which decays to zero over some time scale. 

We follow a similar procedure as in Section [TTJ By cal- 
culating commutation relations, we find; 

id t a z k = 4$^ {°t°J - °t a k) + 4 (% °t - 4°k) 

(29) 

where r] k — | ± iry/). 



idtcrt 



2Bf*+ + 2 V +** k + 2 E J&M (30) 



where Bf = B k + £ . , fe jjtf 
netic field at site k. Also 



idtal 



2Bf al- 2^^1-2^ J fk^] 



Tjfo is the effective mag- 
(31) 



Finally, 

id t (a+a- 



:2ABeff (T + c7 - 



J 



2 V [J^l^-J^crta-] (32) 



where AB|f 



ml 



rf — if,, and j ^ k. Analagous to 



Eqs. (|5|) and $Q of Scction|TTl ABS? consists of a constant 



part, given by <5Bfcj [see Eq. (|SJ)], and a fluctuating part, 
which is now given by 

AB^ ct (t) = (t) + r,i (t) n* (t) , (33) 



compared with Eq. (|6|). We wish to integrate Eqs. (|30[) . 
(f3"TT) . and (|3"2"T) . and thereby find a closed form for the 
time evolution of cr| from Eq. (l2"9l) . 

We start with Eqs. (|3U|) . and (|31|) and apply the same 
logic as in Section [TT] regarding the self-averaging nature 
of the summations (due to a fluctuating Larmor preces- 
sion frequency). What is left can easily be integrated to 
give 



*fc(t) = T2i 



,±2if*Bf (r)dr + 



ds- 



± p ±2if*Bf(r)dT 



c,, e 



(34) 
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where c k = <r^ (t =0) gives the contribution from the 
initial conditions. Looking now at Eq. (|32p . and ignoring 
the summation term, we find 



[Jj- k + 2rg (s)a+(s)] a* k (s)\ds + c jk A kj (0,t) 

(35) 

where A' jk (s,t) = e 21 ^M*", and c jk = afa k (Q) is 
the initial condition. Substituting Eq. (|55|) into Eq. d2H) , 
we find that the terms 2?7+ (s)a k (s) and 2r\ k ~ (s)<jj~(s) 
within the square parentheses of Eq. (|35|) arc summed 
over, and hence can be ignored, due to our self-averaging 
approximation. We then proceed with the same mean- 
field approximation as in Section HH this time replacing 
A'j k (s,t) —> (A'j k (s,t)), which is again assumed to be a 
Gaussian random variable, such that 

(A' k -(s,t)) = e r 2t5Bk ^ t ~ s) e~ A K ' k ' {T){]t ~ s] ~ T)dT 



where 



(36) 

K' kj (t-t') = {AB^(t)AB^(t')) 

= K kj {t-t') + 2k(t-t'). (37) 



Proceeding in this way, Eq. (|29|) for the time evolution 
of <t| becomes, 



&k 



*l {v k (t)e 2i ^ B ^ dT vUs)+ 

4^ 2lSlBf{T)dT ^^)yi^)ds (38) 



where 



ft (t) =4iJ2 Jjh hkA' kj (0, t) - c kj A' jk (0, t)] (39) 



and 
Vk(t) = a 



(40) 

are both noise terms, arising from the initial conditions 
of ato k and <r - respectively. The new rate, T'- k is now 



given by 

?' jk =4(J^) ; 



exp 



2i<5B 



jkS- 



\s\ 



4/ ^.( M )(| S |- M )d M 
'o 



ds, 



(41) 



where we have employed the Markov approximation, to 
remove [<r| (t) — crj: (f )] from the integral, and used the 
quickly decaying property of (A'j k (s, t)) to extend the 
upper and lower limits of the integral to ±oo. 

The integral term in Eq. (|38|l can be greatly simplified 
by replacing the terms in the curly parentheses by their 
average value. This approximation is consistent with an 
assumption of the differing time scales between fluctu- 
ating local magnetic fields at site k, and the individual 
dynamics of a single spin at site k. When a large number 
of individual spins contribute to the local effective field 
B k at site k (as is the case for systems with long range 
interactions or high dimensionality) the fluctuations will 
appear Gaussian, and the term in Eq. (f38|) involving the 
integral, becomes 

-8/ A(t- s)cos[2B k (t- s)] x 
Jo 

e -4/l'-*l G k (T)[\t- S \-T]dT a Z^ ds (42) 



where 



Gk(r) = E E i4«(0)<W) +A(r). 

m^k n^k 



(43) 



The term preceeding a k (s) in Eq. (J32) decays much faster 
than the evolution of crl(s), so we can apply the Markov 
approximation cr k (s) — > cr k (t), and extending the upper 
and lower limits of integration to ±oo wc find 

d t {al) = £ T' jk (a] at)-T k (al) + ft + Vk (44) 



where 



A(s) cos {2B k s) e - 4 ^ S ' G k (n)[\s\-nW ds 



.(45) 

gives a new rate at which the spin direction at site k 
relaxes down into a completely random orientation of ei- 
ther ±1. This relaxation mechanism is entirely due to 
the fluctuating external magnetic field terms; Tjf(t), in 
the Hamiltonian of Eq. 



A. Example: white noise 

If we consider the following simple example 

A(t) = \6(t) (46) 

then wc find = &6 jk \8(t - t') and T fe = 4A. 

We wish to examine two different limiting cases; 

1. y/X < J lk and J lk < B t - B k 

2. J lk < v/A and J ife < Bi - B k 

In case 1. the external noise is sufficiently weak that 
the relaxation time-scale is essentially infinite, in which 
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case we can set (^(i) "™^')) : 
K' jk (t) = n jk + 2X5(t), where 



5 mn . In this way wc find 



m^{j,k} 



Til _ Til 

J mk J mj 



(47) 



Continuing with the calculation, we find the following 
expression for the rate; 



r' 



jk 



8(4) 2 



ds cos [2(B k 



Bj)s] c 



-2[K jfc s 2 +4As] 



(48) 

This integral can be expanded to first order in the small 
parameter, to give 



r' 



jk 



i(Jjk) ds cos [25B ]k s]e- 
Jo 



4(4) 2 



'7r exp 



( (SB ]k ) 2 \ 
\ 2 Kjk J 



4V26B jk \Fv 



4A 

Kjk 
SB ]k 



3/2 



(49) 



where Fn(x) = e~ x Jq e v dy is Dawsons integral^ 2 -. This 
result is shown in the solid lines of Fig. [3] for the case of 
the Heisenberg model on a cubic lattice (as discussed in 
Section [TTTJ) . We can further approximate Dawsons inte- 
gral, in the case of a large gradient 5Bjk H> yJ2Kj k , to 
give Fd(x) ~ + 4^3- +0{x~ 5 ), for large From this 
we find the asymptotic behaviour of the rate 



r' 

L jk 



16(^ fc ) 2 A 



5B% 



(50) 



valid when 5Bj k ^S> Jf k - 

In case 2. the external noise is sufficiently strong, that 
it dominates over the interaction-induced spin-diffusion 
process. We can then approximate Eq. (|4"4"]) as 



9 t (4)~-T fc (a|)+77 fc . 



(51) 



In this case one would observe exponential decay in the 
autocorrelation function (due to the noise term r\k) given 
by 



(52) 



This leads to K' jk (t) = 2X5 (t) + n^e'^ ~ 2X5{t), 
which gives us the following expression for the rate; 



r' = 16(J; 



A 



3k ' SB% 



16A 2 ' 



(53) 



This result is shown in the dashed lines of Fig. [3] for 
the case of the Heisenberg model on a cubic lattice (as 
discussed in Section UTTT) . 
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FIG. 3: Plotting the diffusion rate T'j k /J as a function of 
magnetic field gradient (Bj — B k )/J for a variety of different 
values of A/ J. The solid lines show the result in Eq. (|49p for 
case 1. The dashed lines show the result of Eq. (|53[) for case 2. 
The actual model is taken to be the same as the Heisenberg 
model discussed in Section [XTTl 



Thus, in both cases 1. and 2. we find that the rate 
now decays as the inverse of the gradient squared; ~ 

\Jjj./5Bjkj ■ This provides a huge contrast with the 
noiseless situation of Section UH where the rate decays 



as 



exp 



- [Jlk/SBp 



The presence of the noise 



provides a means for spin diffusion to occur over a much 
faster time-scale (in the presence of a strong external 
magnetic field gradient). 



V. DISCUSSION AND SUMMARY 

In Section |TT] of this article we have derived a dynami- 
cal mean-field theory for systems of spin-half particles on 
a lattice, in the presence of a nonuniform, external mag- 
netic field. The theory is applicable in the case where the 
magnetic field gradient between two lattice sites is large 
compared to the interactions. Additionally, the number 
of interacting pairs should be large (as is the case for sys- 
tems of high dimensionality or long range interactions). 
This condition is necessary to ensure the fluctuations of 
the effective field at each lattice site are Gaussian (the 
central limit theorem) . One of the most notable approxi- 
mations we made in deriving this theory of spin diffusion 
was the exclusion of the summation in Eq. ©. With 
this sum excluded, we were able to derive a solution to 
Eq. §5§ , shown in Eq. §8§ . We can use this expression for 
a^a k (t), to estimate the size of the summation term in 
Eq. ([3]), and thus estimate the error in this approxima- 
tion. 

First, we note from Eqs. (fTTj) and (fTU|) . the size of 
(j^a~f(t) is roughly Tjk/{2J^ k ). Thus, if we substi- 



tute our expression for crfcrf back into Eq. ((3]), we see 
that the size of the summation term is approximately 
max{r„i, T nk } where n runs over lattice sites which are 
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mutual neighbors of sites i and k. Assuming a certain 
level of isotropy exists within the system, we conclude 
that, provided J^. T^, for all interacting pairs i and 
k, the exclusion of the summation in Eq. ([3]) is justified. 
With all conditions satisfied, the equation of motion for 
the z-component of the individual spins is a Langevin 
equation with additive noise, see Eqs. (|12[) . 

If the condition 3> were not satisfied, and the 
summation in Eq. ([3]) could not be justifiably ignored, 
we would expect a similar analysis to be possible. The 
summation term would manifest as multiplicative noise 
in the coefficients of the Langevin equation (fT2j) . as 
well as the additive noise which we have derived. Further 
work on this issue however, is still in progress, and the 
details deferred to a future publication. 

The model can be described in terms of simple phys- 
ical principles, as illustrated in Fig. [T] Interactions 
between sites i and j can cause spin flip-flopping, i.e. 
\if3i.) luit)- This process occurs when sites i and j 
have opposite spin, and does not conserve energy when 
the external field gradient is nonzero (due to the different 
Zccman energies). These sites i and j however, also inter- 
act with all other neighboring lattice sites (the number of 
which is assumed to be large). A crucial approximation 
in our model is to treat all remaining sites as composing 
an effective bath, or rapidly fluctuating environment in 
which sites i and j inhabit [see Fig. [T](b)]. In this way, 
one can derive the rate at which the spin flip-flopping 
occurs (we have labelled this quantity Fy), and natu- 
rally it will depend on the bath parameters. To be more 
specific, it depends on the correlation functions between 
neighboring sites within the bath. The final step then is 
to determine the rate that is self-consistent with the 
bath, i.e. the value of Fy which yields the same corre- 
lation function between neighboring sites, as that from 
which it was derived. 

We find the rate decays very quickly with increas- 
ing field gradient. Equation (fl~6|) predicts the rate decays 
in the same way as a Gaussian distribution. From a nu- 
merical study of the cubic Heiscnbcrg lattice (presented 
in Section llllj) . we expect this prediction to be accurate 
for B, t - B 3 > lOJy (see Fig. [2]). This result implies that 
the observation of spin diffusion in systems with a very 
strong magnetic field gradient is likely to be difficult as 
the diffusion time-scales would be very large. 

However, in Section [IV] we studied the influence of ex- 
ternal noise on this rate. The presence of the external 
noise turns out to be favourable for increasing the diffu- 
sion rates. We made use of the same set of assumptions 
in deriving a second Langevin equation [see Eq. (|44[) ]. 
In contrast to Section HH the Langevin equation now in- 
cludes a decay-constant, denoted which relaxes the 
system down into a state where the orientation of the 
magnetic moment is completely random, i.e. (cr|) = 0. 
Spin flip-flops still occur in the system, and the rate at 
which they occur; T^-, is affected by the noise. As a gen- 
eral rule, the rate T'^ increases with increasing noise, as is 
illustrated in Fig.|3J In the limiting case where the exter- 



nal noise is far greater than both the interaction coupling 
and the external field gradient, we find the rate de- 
cays in the same way as a Cauchy-Lorcntz distribution, 
see Eq. (|53p . This predicted increase in the rate may help 
to explain experiments where diffusion has purportedly 
been observed in systems with very large magnetic field 
gradients. 
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Appendix A: Perturbation theory for the two-body 
problem in a fluctuating external field 

Consider the following time dependent Hamiltonian 
describing two spin-half particles located at sites 1 and 
2, interacting via an exchange interaction, 

H = B 1 (t)c7 z 1 +B 2 (t)al + J(afa^ + a\al + a{a z 2 ) (Al) 

where the external magnetic field Bi(t) = Bi + bi(t) con- 
sists of a constant part and a fluctuating part. We wish 
to calculate the probability of the spins flip-flopping in 
time t, that is 



p(t)= (l 4 2 t |E/(t,0)|l t 2 4 ) 



(A2) 



where U(t, 0) is the time evolution operator for H . Split- 
ting the full Hamiltonian up into a noninteracting and 
an interacting part, 



H 
V 



B 2 (t)a z 



'l u 2 



o\o z 2 ) 



(A3) 
(A4) 



and moving to the interaction picture; |\f's(i)) 
U {t,0)\^i(t)) where U (t,0) = cxp -i f* H [r)d 
is the time evolution operator of the noninteracting 
Hamiltonian. Defining V/(i) = Uo(t,to)VUo(to,t) and 

Ui(t, 0) cxp + —i J Q Vi(r)dr , where cxp + denotes the 

usual time-ordered Dyson series^ (appropriate for non- 
commuting [Vi(ti),Vi(t 2 )] ^ when t\ ^ t 2 ). 

Using this standard formalism, wc approximate the full 
time evolution operator as, 



U(t,0) = U o (t,0)Ui(t,0) 



U o (t,0) 



l-i I Vi{r)dT 
o 



(A5) 
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by truncating the Dyson series for Uj. Substituting this 
approximation into Eq. (|A2[) and working through the 
calculation in a straight-forward manner we arrive at 



p(f) = 4 J 2 



dry [ 
o Jo 



dr e -2i(B 1 -B 2 )(T 1 -T 2 )+2if;2 8b(s)ds 



where Sb(s) = bi(s) — &2(s)- At this point it is convenient 
to take an average over the fluctuating component of the 

i c ii 2i r 2 Sb(s)ds , 2i P 2 5b(s)ds\ . 

external held, e Jt i — > (e Jt i ). Assuming 

these fluctuations are Gaussian, and time-translationally 
invariant, we hnd 



dT2e -2,:(B 1 -S 2 )(r 1 -r 2 ) x 



p(t) =4 J 2 / dn / 
Jo Jo 

p -i f^ T1 - T21 ds(Sb(Q)Sb(s))(.\r 1 -T 2 \~s) 



4J 



dr (*— M), 



-2i(Bi-B 2 )r, 



D -4/ r J T| d s (56(0)5b( S ))(|r|- s ) 



(A6) 



In the limit then, where the time t is much larger 
than the time scale over which the final term 
e - 4 /d T ' ds(Sb(o)Sb(s))(\r\-s) decayS; the probability be- 
comes 



p(t) =4tJ 2 



exp 

1*1 



2i(B 1 -B 2 )s~ 
(Sb(0)Sb(s)) (|s| - ft) dfj, 



This probability in Eq. (|A7[) should be compared to 
the rate at which spin flips are predicted to occur from 
Eq. (|I3[) in Section UT] In making this comparison, we see 
that the approximations we have applied in deriving the 
equation of motion (fT2"]) for cr| amount to treating all sites 
other j and k as composing an effective bath (equivalent 
to a fluctuating external field). 



Appendix B: Numerical algorithm for solving the 
integral equation 

For our particularly simple choice of B(r) = -^={r x + 

r y + r z), the integral equation we must solve is simply 
Eq. ([24]) with K(fi) given by Eq. ([26]). From this equa- 
tion, we wish to determine T as a function of bo, the de- 
pendence on J can be removed, by switching to variables 
T = T/J and b = bo/ J, such that we have 

d se ^ e -8/r^e- 6r ^/(2f>)(| S |-^ 

-oo 

We then search for a root of this equation, by iterating 

f(n+l) =4 / >00 dse ^7§ f e -8/l 3| rfMe- 6f( "^/(2f<") M )(| s |- Al )^ 

(Bl) 



for n ~ 0, 1, 2, . . . up to convergence, which in our case 
was chosen to be |f( n+1 ) - f»| < 1CT 4 . In order to 
choose a reasonable initial prediction for we begin 
the algorithm at bo = 20, and define 



f(o) 



20 



(B2) 



Once the algorithm has converged, we decrease bo by a 
small amount and use our previous prediction for T as 
our new T^°\ 



Appendix C: The issue regarding 
convergence/divergence of T as 6o -> 

As the field gradient decreases in a particular direction, 
the rate at which spin flip-flops occur in that particular 
direction increases, see Figure [5] It is not clear, a priori, 
that the rate will remain finite in the limit of vanishing 
gradient. Consider, for example, the RHS of Eq. ([24]) 
(and set J = 1). We can rewrite this in terms of the 
Fourier transform of K(fi) = j- J e~ lu) ^K{oj)doj, and we 
are only interested in the case where bo = 0, so we find, 



RHS(r) = 4 / exp 



ds (A7) Next we define 7 = fiT, in which case 



ds. 
(CI) 



dj 



-e r e 



-by 



/(27) 



r vr 



(C2) 



where H(x) = J R dj e^e- 6 ^ f{2j). The RHS therefore 
becomes, 



RHS(r) = 4 f exp -- f dr}H(r]) 



1 - cos(r?7|s|) 



rv 



ds. 
"(C3) 

where we defined 77 = uj/T. Now we make the assumption 
that r does become very large, in this limit we find 



I-cos(r?7|5|) ir\s\ 

rv > ~ 6{v) 



and therefore 



RHS(r) 



2r 



(C4) 



(C5) 



for large T. The quantity H(0) can be calculated nu- 
merically to be H(0) ~ 2.33, thereby indicating that the 
slope of the RHS is < 1 for large T. From this we con- 
clude that a finite value of T will exist at 60 = which 
satisfies Eq. ([2"3]) . 
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